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Abstract 

Randomness and frustration are considered to be the key ingredients for the existence of spin 
glass (SG) phase. In a canonical system, these ingredients are realized by the random mixture of 
ferromagnetic (FM) and antiferromagnetic (AF) couplings. The study by Bartolozzi et al. [Phys. 
Rev. B73, 224419 (2006)] who observed the presence of SG phase on the AF Ising model on scale 
free network (SFN) is stimulating. It is a new type of SG system where randomness and frustration 
are not caused by the presence of FM and AF couplings. To further elaborate this type of system, 
here we study Heisenberg model on AF SFN and search for the SG phase. The canonical SG 
Heisenberg model is not observed in d-dimensional regular lattices for [d < 3). We can make an 
analogy for the connectivity density (m) of SFN with the dimensionality of the regular lattice. It 
should be plausible to find the critical value of m for the existence of SG behaviour, analogous to the 
lower critical dimension {di) for the canonical SG systems. Here we study system with m = 2,3,4 
and 5. We used Replica Exchange algorithm of Monte Carlo Method and calculated the SG order 
parameter. We observed SG phase for each value of m and estimated its corersponding critical 
temperature. 
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I. INTRODUCTION 


Spin glass (SG) is one of the most complex systems in condensed matter physics and has 
been intensively studied in the last four decades[l-6|. It is a randomly frustrated magnetic 
system with frozen disordered spin orientation at low temperatures. This unusual conhgu- 
ration is regarded as a temporal ordered phase j^, different from the spatially ordered phase 
found in regular magnets. The complexity of the system is due to the presence of frustration 
and randomness, which are the key ingredients for the existence of SG phase. Frustration is 
a state where spins can not hnd hxed orientations to fully satisfy all the interactions with 
their neighboring spins. This can be caused either by the conflicting interaction between FM 
and AF couplings, or between among AF coupling due to the topological factors. Frustra¬ 
tion alone can not lead a system to an SG phase, hrmly exemplihed by the fully frustrated 


AF planar spin systems which have spatially ordered phase at low temperatures 


,l9|. 


Most SGs studied are canonical system where both FM and AF couplings exist. The 

n n 

examples of these are Sherrington-Kirkpatrick model[^, Edward-Anderson model[3[ and p- 


spin interaction model[^. Bartolozzi et al. hrst reported SG behavior of the Ising model 
with AF interaction on scale free network (SFN)[10|. This is a new type of SG system 
without random distribution of FM and AF couplings. The nodes of the network do not 
have homogeneous number of links and frustration is fully due to the topological factor. 
The work has brought a new insight into the study of SGs, suggesting that the irregular 
connectivity can also be one of the ingredients of SGs, different from the previous notion 
insisting the presence of random mixture of FM and AF interactions. 

An SFN consists of abundance of triangular units on which spins are frustrated if the 
couplings are AF. While the frustration is caused by the topological factor, the randomness 
is due to irregular connectivity. There are some vertices having very large number of connec¬ 
tions, acting as the hubs as in the internet connection. In fact, the structure of the internet 
follows a scale free behavior. The AF Ising model on random networks without a scale free 

n 

behavior was also reported to exhibit SG phase [1I[. We believe that random connectivity 
can generally be an alternative for the usual randomness, together with frustration, as the 
ingredients of SGs. 

Here we study the Heisenberg SG model on SFN. The prevalent controversy on existence 
of SG phase for the 3D canonical Heisenberg model ^ over the last three decades is one of the 












main motivations. Most SGs studied, such as CuMn, AgMn or CdMnTe, are Heisenberg-like 
systems. Early numerical study by Coluzzi observed SG phase transition on 4D Heisenberg 
model. Kawamura and Nishikawa pointed out the absence of Heisenberg SG phase on D- 
dimensional space, for D < 3 j^. These works suggested that the lower critical dimension 
di might be a fractional number, 3 < d/ < 4. 

For a regular lattice, spatial dimension is related to the coordination number, i.e., the 
number of neighbors of each spin. We can associate the coordination number with con¬ 
nectivity density m, which is the average number of links, of SFNs. Due to the presence 
of short-cut between spins, the notion of spatial dimension is not an appropriate term for 
SFNs. Nonetheless, an SFN can still be associated with a high dimensional regular lattice. 
Therefore, we assume that there may also be a critical value of m, analogous to di of regular 
lattice, for the random and SFNs. The paper is organized as follows: Section HIl describes 
the models and method. The results are discussed in Section m Section nYi is devouted for 
the summary and concluding remarks. 


II. MODEL AND METHOD OF SIMULATION 


The Heisenberg model on an SFN can be written with the following Hamiltonian, 


h = -jY.' 


( 1 ) 


(u> 


where for an AF system, the coupling constant is set to be negative (J < 0); and s* are 
the Heisenberg spins residing on the nodes of the network. The summation is performed 
over all directly connected neighbors. In an SFN, the number of neighbors of each spin is 
not homogeneous. It is distinguished from a random network due to its scale free behavior, 
P{k) = k~'^, where k is the number of links of each node and 7 the decay exponent of its 


link distribution 


12| . Networks with large 7 have very famous nodes (or hubs), i.e., those 


having direct links to most other nodes. This type of network found many realizations in 
real world, from World Wide Webs, power grids, neural and cellular networks, till routers 
of the internet and citation network of scientists. 

In addition, an SFN is also characterized by a clustering coefficient, C, defined as the 
average of local clustering, Q, 

C,= 




Zi{Zi - 1 ) 


( 2 ) 




where Zi is the total number of nodes linked to the z-th site and yi is the total number of 
links connecting those nodes. Both Q and C he in the interval [0,1]. For a fully connected 
network, C = Cj = 1, all nodes connect to each other. The parameter C is related to the 
number or the density of triangles. Because spins on the vertices of triangular units are 
frustrated, C is also related to the frustration degree of the network. More comprehensive 
review on the ubiquitous of SFNs can be found for example in 13]. 

Thermal averages of the physical quantities of interest are calculated using the Replica 


Exchange algorithm of MC method 


16| . This algorithm is implemented to overcome the 


slow dynamics due to the presence of local minima in the energy landscape of the system. 
Slow dynamics is a phenomenon commonly found in dealing with SGs, where a random 
walker can easily get trapped in one particular local minimum. It is an extended Metropolis 
algorithm where a system is duplicated into K replicas. All replicas are simulated in parallel 
and each is in equilibrium with a heat bath of an inverse temperature. Given K inverse 
temperatures, • • • ,[3k, the probability distribution of finding the whole system in a 

state {X} = {Xi, X 2 ,..., Xk} is given by, 

K 

P{{X,/3})=l[P{X^,^,n), (3) 

m=l 

with 

P{X^,M = Z{/3J-^exp{-P^H{Xj), (4) 

and Z{j3m) is the partition function of the m-th replica. We can define an exchange matrix 
between replicas as hF(Xm,/dn), which is the probability of switching configuration 
Xm at the temperature /3m with configuration X„ at /3„. 

In order to keep the entire system at equilibrium, by using the detailed balance condition 


P(. . . , {Xm, ^m}, . . . , {X„, /34, . . .) • VF(Xm, /3m|X,„ (3^) 
= F(. . . , {X„, /3m}, . . . , {Xm, /3J, . . .) • W(X^, /3m|Xm, /3„), 


(5) 


along with Eq. (|1]), we have 


W(X^,/3jX^,/3„ 

W(X^,/3jX^,/3„ 


= exp(-A), 


( 6 ) 


where A = {H{Xm) — H{Xn))[j3n — (3m)- With the above constraint we can choose the 
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FIG. 1: Energy time series of system with L = 4096 and m = 5 at T = 0.2 and 0.5. 
matrix coefficients according to the standard Metropolis method, therefore 

f 1 if A < 0, 

W{X^,(5m\Xn,l^n)= { ~ (7) 

I exp(—A) if Otherwise. 

As the acceptance ratio decays exponentially with (/3„ — jSm), the exchange is performed 
only to the replicas next to each other, i.e., W{Xm, l3m\Xm+i) The replica exchange 

method is extremely efficient for simulating systems such as SGs. This method has been 
widely implemented in many complex systems, including the AF-SFN Ising Sojl^. In the 
next section, we discuss the search for SG behaviour of Heisenberg antiferrmagnetic system 
on SFN. 

III. RESULTS AND DISCUSSION 


A. Equilibrium Behavior 


We have simulated AF Heisenberg model on SFN with various connectivity densities, 
m = 2,3,4 and 5. For each density, we consider several system sizes, L = 1024,2048, 
and 4096, which are the number of spins. Since this is a random system, we took many 
realizations of the networks for each system size, then averaged the results over the number 
of realizations. Each realization corresponds to one particular connectivity distribution. We 
have to take reasonable number of realizations Ah for the better statistics of the results. 


Previous study on Ising system took = 1000 reaIizations[l^. Here, due to the less 
fluctuation of the results from different realizations, we took moderate number of realizations, 
i.e.. Nr = 500. 

































To check for the equilibrium behavior of the systems we evaluated the energy time series. 
In Monte Carlo (MC) simulation, time is assigned as a series of MC steps (MCSs). One 
MCS is dehned as visiting each spin once, either randomly or consecutively, and perform¬ 
ing a prescribed spin npdate, i.e.. Metropolis update. We performed M MCSs for each 
temperature and took N samples ont of M MCSs. The time series plot of energy for two 
different temperatnres for system size L = 4096 is shown in Fig. [H The phenomenon of 
slow dynamics at lower temperatnre (T = 0.2) is shown in Fig. [1] where the average valne 
for the hrst half of the plot is larger than that of the rest half. To overcome this problem, 
the eqnilibration process was increased np to 8000 MCSs. 




FIG. 2 : Temperature dependence of (a) specific heat and (b) SG order parameter for different sizes 


with m =2,3,4 and 5. The solid lines are guides to the eye for the clarity of peaks. 


We also calculated the specific heat which is defined as follows 

a = ^ ({E^) - {Ef) (8) 

where L, k and {E) are respectively the nnmber of nodes, Boltzmann constant and the 
ensemble average of energy. The temperatnre dependence of the specihc heat for various 
system sizes is shown in Fig. [2](a), where the statistical errors are comparable to the size of 
the symbols. Although the specihc heat plot exhibits no singnlarity, there is a maximum 
value at intermediate temperatnres. As indicated, the peaks of the specihc heat shift to 
higher temperatures as the increase in m. The presence of peak may signify the existence of 
phase transition. A clear sign of the SG phase can be observed from the plot of overlapping 
order parameter which will be presented in the next subsection. 








FIG. 3: Binder parameter of qEA for various values of m. Solid lines are guide to the eyes. 

B. Spin Glass Order Parameter 


To search for the SG phase transition, we calcnlate SG order parameter dehned as follows 


qEA 


i 


(9) 


where s* is the spin on node f-th while a and /3 denote two sets of replicas. The origin of this 
qnantity is from the scalar prodnct of the vector spins. A scalar product of two vectors will 
be maximum (zero) if they are parallel (perpendicular). This idea is implemented to capture 
the frozeness of the spin conhguration. If a system is frozen, the spin conhgurations of two 
replicas with the same inverse temperature from different sets will be more or less similar, 
regardless of globally rotational invariant. As explained in Sec. [TTl system is replicated into 
K replicas, each belongs to inverse temperature. For the sake of qea calculation, the whole 
replicas are duplicated, therefore we have two sets of replicas. Each set contains K replicas 
which are exchanged during the simulation. Only replicas from the same set were exchanged 
during the simulation. 

The overlapping parameter will give hnite value if system in SG phase. This applies 
for any SG model, including the Ising and the Heisenberg model. For Ising model, qea 
is simply the multiplication of the overlapped spins. In contrast for the Heisenberg model 
where spins are allowed to rotate in any direction, we take the tensor product instead of the 
scalar product, resulting qea with nine components. The plot of temperature dependence 
of qea for various different system sizes is shown in Fig. [2](b). As indicated, this parameter 
increases as temperature decreases, which is the evidence for the existence of SG phase at 
lower temperatures. 















To clarify that this a true SG phase, we calculate the cumulant ratio (Binder parameter) 
of qea defined as follows 


Ur = 


11 - 


9(g^) 

(g2)2 


( 10 ) 


The plot oi Ul for different system sizes is shown in Fig. |3l It is clear that there is a 
single crossing point, emphasizing the existence of SG phase transition. However, due to the 
anomaly of the crossing patterm which is different from the standard ones, such as in the 
previous study on Ising model [l0|, we did not perform scaling plot oiUi for the estimate of 
critical temperature and exponent. Instead, we roughly estimate the critical temperature for 
each connectivity density as listed in Table 1, where numbers in bracket are the uncertainty 
for the last digits. There is a systematic increase of critical temperature as the connectivity 
density increases. This is consistent with the Mean Field theory argument, where an ordered 
phase for a system with large connectivity tend to be more stable. 


IV. SUMMARY AND CONCLUSION 

In summary, we have studied AF Heisenberg model on scale free network, using Replica 
Exchange MG method. We simulated several different connectivity densities (m = 2,3,4 
and 5) and calculated such physical quantities as ensemble average of energy, the specific 
heat, the overlapping parameter and its cumulant ratio (Binder parameter). A sign for finite 
SG phase transition was observed for all values of m. There is a systematic increase of Tc 
as m increases. This is related to the fact that systems with large connectivity tend to be 
more robust against thermal fluctuation. The existence of SG phase even in the system with 
m = 2 is generally in a good agreement with the case for canonical systems where Heisenberg 
SG phase was clearly observed in system with large dimension, e.g. the 4D system[M|. It 
is therefore interesting to search for the lower critical value of m by studying systems with 
fractional value of m. We will consider such system in our future study. 
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TABLE I; The estimate of critical temperature for each value of m. 


Connectivity density 

Tc 

(m) 


5 

0.44(3) 

4 

0.35(5) 

3 

0.24(6) 

2 

0.19(6) 
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